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A METHOD OF ESTIMATING PRECIPITATION CHARACTERISTICS 



[0001] The present invention relates to a method of estimating precipitation characteristics 
and in particular precipitation rate for solid precipitation. 

[0002] It is known from the state of the art that the characteristics of a radar image can be 
used to estimate rainfall characteristics. In particular, European Patent EP 1 049 944 describes a 
technique for estimating rainfall using a radar. In that patent, provision is made for the following 
steps: 

using that dual-polarized radar to measure the differential phase shift (Odp) and 
the apparent reflectivity Z in at least one of the H (horizontal) or V (vertical) polariza- 
tions and over a given range [ro, ri] of path radii relative to said radar; 

determining an estimate of the value No* representative of the size distribution of 
the rain drops, on the basis of the differential phase shift difference over the range ro to 
ri, and on the basis of an integral of a function of the apparent reflectivity Z, over the 
range [ro, ri]; and 

deducing the value of the precipitation rate at a point on the basis of No* and of 
the apparent reflectivity at said point. 
[0003] A method of estimating rainfall is also known from PCT Patent WO03007016 which 
describes a method of estimating a precipitation rate by means of a dual-polarized radar, that 
method being characterized by the following steps: 

using that dual-polarized radar to measure the differential phase shift Odp and the 
attenuated reflectivity Z in at least one of the H or V polarizations and over a given range 
[rl, rO] of path radii r relative to said radar; 



determining an estimate of the value K(rO) of the specific attenuation at rO on the 
basis of the profile of the attenuated reflectivity measured in this way, and on the basis of 
the differential phase shift difference over the range rO to rl; determining an estimate 
K(r) of the specific attenuation at r as a function of the attenuation K(rO) determined in 
this way and of the attenuated reflectivity profile Z(r); and determining the precipitation 
rate R(r) once K(r) is known. 
[0004] Those various solutions make it possible to characterize liquid precipitation and to 
estimate precisely rainfall rate (in millimeters per hour (mm/h), but they do not make it possible 
to estimate the characteristics of solid precipitation such as snow. 

[0005] An object of the invention is to provide a novel method that remedies that drawback. 
[0006] To this end, and in its most general acceptation, the invention provides a method of 
estimating precipitation characteristics and in particular precipitation rate for solid precipitation, 
said method comprising an acquisition and processing step consisting in acquiring a radar image 
including at least a vertical plane of a precipitation zone, and in processing a vertical profile so as 
to deliver digital signals representative of reflectivity in the vertical direction h, said method be- 
ing characterized in that it further comprises an integration step consisting in integrating said sig- 
nals representative of reflectivity so as to deliver a signal representative of the profile in the ver- 
tical plane of a the mean particle diameter weighted by the mass of each particle, and a deter- 
mination step consisting in determining the concentration of the solid particles on the basis of the 
signals computed in the preceding steps. 

[0007] Preferably, the integration step consists in determining the variable Z(h) of the radar 
observable in mm 6 /m 3 as a function of the altitude h on the basis of said radar image, and in de- 
termining said mean diameter Dm(h) of the particles by resolving the following equation: 
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where: 



Z is the radar observable to be inverted 
D m is in meters (m); 

a and b are coefficients specific to particles of the "aggregate" type; for example, the co- 
efficient a is equal to 35184 and the coefficient b is equal to 3.16; 

keg-is the coefficient of effectiveness of the aggregation process to be adjusted, said coef- 
ficient keff being equal to 0.3. 

[0008] Integration of (2) requires an integration boundary condition. Advantageously, said 
integration boundary condition is determined so that the value D m (h) at the top of the cloud cor- 
responds to the predetermined value for the total number of particles at the top of the cloud. 
[0009] In an advantageous implementation, the profile of the total number of particles n t (h) 
is determined by the following equation: 

n T (h)=x.Z(h)/D m (h) 6 
where x is equal to 25.4 10~ 18 . 

[0010] In another implementation, the meteorological parameter No(h) is determined by the 
following equation: 



where y is equal to 1 02 1 0" . 

[0011] In a third implementation, the meteorological parameter corresponding to the profile 
of the ice water content IWC(h) (in g/m 3 ] is determined by the following equation: 



N 0 (h)=y.Z(h)/D m (h) 7 



IWC(h)=wZ(h)/D m (h) 3 



where w is equal to 1.25 10 
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[0012] In a fourth implementation, the meteorological parameter corresponding to the profile 
of the solid precipitation rate R(h) (mm/h equivalent melted) is determined by the following 
equation: 

R(h)=r.Z(h)/D m (h) 235 
where r is equal to 4.698 10~ 10 

[0013] The invention will be better understood on reading the following description, given 
with reference to a non-limiting implementation. 

- Figure 1 shows an example of a vertical profile of Z to be inverted (in this example, the 
isotherm 0°C is at ground level); 

- Figure 2 shows the profile of Dm resulting from Z being inverted by the aggregation 
model compared with the conventional estimator; 

- Figure 3 shows the profiles of NO and of nT resulting from Z being inverted by the ag- 
gregation model compared with the conventional hypothesis and observations; 

- Figure 4 shows the vertical profile of Z to be inverted; 

- Figure 5 shows the profile of Dm resulting from Z being inverted by the aggregation 
model compared with the conventional estimator; 

- Figure 6 shows the profiles of NO and nT resulting from the Z being inverted by the ag- 
gregation model compared with the conventional hypothesis and observations; 

- Figure 7 shows the profile of R resulting from Z being inverted by the aggregation 
model compared with the conventional estimator; 

- Figure 8 shows the sensitivity of the retrieval of the equivalent precipitation rate (mm/h) 
to the particle density law p(D) oc D-y; and 
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- Figure 9 shows the sensitivity of the retrieval of the vertical Doppler velocity to the par- 
ticle density law p(D) <x D-y. 
[0014] The method applies essentially to stratiform precipitation. It considers that the ice- 
forming nuclei are activated only at highly negative temperatures, i.e. at the top of the cloud. 
The ice crystals formed at high altitude settle and grow as they fall, either by sublimation of the 
ambient saturating water vapor, or by collection and freezing of supermelted cloud water drop- 
lets, or else by random aggregation of their collisions with other ice crystals. Of the three growth 
processes, only aggregation changes the ice particle concentration. The characterization method 
is based essentially on a simplified description of the aggregation mechanism. The invention is 
based on a class of "profiler" methods in that it inverts the vertical profile of reflectivity mea- 
sured in the ice, so as to extract therefrom the vertical profile of the solid precipitation rate. 
[0015] The steps of the inversion method are as follows: 

1 - Particle size distribution expressed in "equivalent melted diameter" is assumed to be 
exponential, i.e.: 

N(D)=N 0 exp(-4D/D m ) (1) 
where N(D) is the concentration of particles per cubic meter (m 3 ) and per diameter range, 
and No and D m are the two parameters that characterize distribution. 

2 - The top h max and the base h min of the layer of solid precipitation are determined; 

a. h max is the maximum altitude of the measured reflectivity profile Z(h). 

b. hmin is either the altitude of the isotherm 0°C if the ground-level temperature is 
positive, or it is ground level if the ground-level temperature is negative. 

3 - The profile of the parameter Dm(h) in the range h max to h min is then determined by re- 
solving the following differential equation: 
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where: 

- Z is the radar observable to be inverted in mmV 3 : 

- D m is in meters (m); 

- a and b are coefficients specific to particles of the "aggregate" type, equal re- 
spectively to 35184 and to 3.16 on the basis of the observations of Locatelli and 
Hobbs(1974); 

- k e ff is the coefficient of effectiveness of the aggregation process to be adjusted 
(the value k e jf= 0.3 seems correct). 

4 - The integration of (2) takes place from the top, where the boundary condition is ex- 
pressed by fixing the total number of particles nr (or the number of ice-forming nuclei 
activated at the top of the cloud). It is possible to take n T (h ma x) = 10 6 m~ 3 , which makes it 
possible to express the boundary condition D m (h max ) as: 

Dm(H max )=25.4 10- l8 (ZQx m Jm(h max )) m (3) 

5 - Once the profile D m (h) from h max to h min has been determined, the profiles of the other 
parameters of interest are computed by the following expressions: 

a. Profile of N 0 : N 0 (h)=102.10' l2 Z(h)/D m (h) 7 

b. Profile of the total number of particles n T (h) [in m' 3 ]: N T (h)=102.10' l2 Z(h) 
/D m (h) 6 

c. Profile of the ice water content IWC(h) (in g/m 3 ]: 

IWC(h)=l .25. 10" 12 Z(h)/D m (h)3 
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d. Profile of the solid precipitation rate R(h) (mm/h equivalent melted). By using 
the terminal fall velocity determined by Locatelli and Hobbs for aggregates: 
[vt=107.6 D 0 65 (D in m)], R(h) is expressed by: 
R(h)=4.698.10" 10 Z(h)/D m (h) 2 35 

[0016] The following description relates to another version of the method of processing for 

determining the precipitation rate. 

[0017] For rainfall, the algorithm used is the ZPHI algorithm which is the subject of French 
Patent Documents FR9800714 and FRO 109206. 

[0018] For solid precipitation, the novel algorithm that has just been developed is the subject 
of the following description. 

[0019] This algorithm for estimating the precipitation rate for solid precipitation belongs to 
the class of "profiler" algorithms in that it inverts the vertical profile of reflectivity measured in 
the ice, so as to extract therefrom the vertical profile of the solid precipitation rate. 
[0020] The algorithm applies essentially to stratiform precipitation. It considers that the ice- 
forming nuclei are activated only at highly negative temperatures, i.e. at the top of the cloud. 
The ice crystals formed at high altitude settle and grow as they fall, either by sublimation of the 
ambient saturating water vapor, or by collection and freezing of supermelted cloud water drop- 
lets, or else by random aggregation of their collisions with other ice crystals. Of the three growth 
processes, only aggregation changes the ice particle concentration. The algorithm is based 
essentially on a simplified description of the aggregation mechanism. The steps of the inversion 
method are as follows: 

Particle size distribution expressed in "equivalent melted diameter" is assumed to 
be exponential, i.e.: 
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N(D)=Noexp(-4D/D m ) (1) 
where D is the equivalent melted diameter of the ice particle; 
N(D) is the concentration of particles per m and per diameter range; and 
N 0 and D m are the two parameters that characterize distribution. 

The top h max and the base h m j n of the layer of solid precipitation are determined; 
hmax is the maximum altitude of the measured reflectivity profile Z(h). 
h m i n is either the altitude of the isotherm 0°C if the ground-level temperature is 
positive, or it is ground level if the ground-level temperature is negative. 
[0021] The profile of the parameter D m (h) in the range h max to h min is then determined by re- 
solving the following differential equation, whose solution can be determined analytically: 

5D m 4 d ke ff J(b,d)a b -(a-2) 1 dZ 

= -24 ZD + D m (2) 

Sh T(4+d) a(a-l) m (a-l)Z dh 

where: 

foo fx r i 

J(b,d)= I exp(-4x)x b I I x d -x ,d |exp(-4x')dx'dx 

Jo JO L J 

describes the collision frequency; 

Z is the radar observable to be inverted in mm 6 m' 3 ; 
D m is in m; 

- a, b, c, and d are coefficients dependent on the density law F(D) for ice particles, 
which law is assumed to vary at D y (where 7 can take values in the range 0.25 to 1.1 
depending on the type of particles). These coefficients are drawn from Mitchell's 
Theory (Journal of Atmospheric Sciences, 53, 12, 1996) for representing the follow- 
ing power laws: 
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v t (D) = cD d (where v, is the terminal fall velocity of the particle of equivalent 
melted diameter D): 

A = aD b (where A is the effective cross-section of the particle of equiva- 
lent melted diameter D); 
With coefficients depending on frequency and on a function of (D), representing 
by a power law the relationship between reflectivity (measured parameter), the concentra- 
tion of the particles, and their mean equivalent diameter: 
Z = aN 0 D m d 

keff is the coefficient of effectiveness of the aggregation process to be adjusted 
(the value 1^= 0.3 seems correct immediately above the isotherm 0°C). 
[0022] The integration of (2) takes place from the top, where the boundary condition is ex- 
pressed by fixing the total number of particles n T (or the number of ice-forming nuclei activated 
at the top of the cloud). It is possible to take n T (h max ) = 10 6 m" 3 , which makes it possible to ex- 
press the boundary condition D m (h max ) as: 

1 

|~i Z(h max ) 1 d-i 

D m (h max )= I 1 (3) 

[0023] Once the profile D m (h) from h max to h min has been determined, the profiles of the other 
parameters of interest are computed by the following expressions: 
Profile of N 0 : 

Z(h) 

N 0 (h)= 

dD m (hf 

Profile of the total number of particles nj(h) [in m" 3 ]: nr(h)=025No(h)D m (h) 
Profile of the ice water content IWC(h) (in g/m 3 ]: 



IWC(h)=(\2212)\tfN 0 (h)D m (hf 
Profile of the solid precipitation rate R(h) (mm/h equivalent melted): 

T(4+d) 

R(h)=(\ .885)10 6 c N 0 (h)D m (h) (4+d) 

4 (4+d) 

Vertical profile of the Doppler velocity V D (h) (in meters per second (m/s)), given 
by a power law dependent on density _(D): 
V D (h)=pD m (h) q 

[0024] Because of the dependency of the coefficients relative to the density law F(D), the re- 
trieval of the precipitation rate depends critically on the parameter 7, as shown by Figure 7. The 
density law thus constitutes a key parameter. 

[0025] By operating in vertical firing, it is possible to measure the Doppler velocity Vo(h), 
and to compare it with V D (h) computed by the algorithm. The Doppler velocity profile is dis- 
criminating, as shown by Figure 8. By successive adjustments, it is thus possible to determine 
the parameter ;yof the law F(D) for which consistency between V D (h) as computed and V D (h) as 
measured is achieved. Figure 4 shows an example of a profile inverted by the algorithm (in this 
example, the isotherm 0°C is at ground level). Figure 5 shows the curves resulting from appli- 
cation of the aggregation model and of the conventional model, making it possible for the profile 
of D M resulting from the inversion of Z by the aggregation model to be compared with the con- 
ventional estimator. 

[0026] Figure 6 shows the curves of the profiles as a function of altitude (m) above the iso- 
therm 0°C. It enables the profiles of N 0 and n T resulting from Z being inverted by the aggre- 
gation model to be compared with the conventional hypothesis and observations. Figure 7 shows 
the profile of R resulting from Z being inverted by the aggregation model compared with the 
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conventional estimator. Figure 8 shows the sensitivity of the retrieval of the equivalent partici- 
pation rate (mm/h) to the particle density law p(D) oc D" Y . Figure 9 shows the sensitivity of the 
retrieval of the vertical Doppler velocity to the particle density law p(D) oc D~ Y . 
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